/*********** NOTES HEADER *********** 

DESCRIPTION: This program estimates the leads and lag merger estimates from the logit
	specification and the marginal effects relative to the period of the consummated merger

Datasets used:
 (1) MatchedSample_Full

Datasets created:
 (1) LL_diab_zip
 (2) LL_hyp_zip 
 (3) LL_death_zip
 
************************************/ 

set more off
set matsize 10000
capture log close
clear all
 
global origData "N:\MedicareClaims-P045601-BE"
global dataIn "N:\MedicareClaims-P045601-BE\Work\hosp_retro\health_out\Data-In\"
global dataOut "N:\MedicareClaims-P045601-BE\Work\hosp_retro\health_out\Data-Out\"
global dpath "N:\MedicareClaims-P045601-BE\Work\ay_data"
global dataProp "N:\MedicareClaims-P045601-BE\Work\hosp_retro\health_out\Data-Out\PropScore\Patients"
global logs "N:\MedicareClaims-P045601-BE\Work\hosp_retro\health_out\Logs\PropScore\Patients"
global dpath "N:\MedicareClaims-P045601-BE\Work\ay_data"
global skapath "N:\MedicareClaims-P045601-BE\Work\ska"

adopath +  "N:/SIL-Common/estout"
adopath +  "N:/SIL-Common/outreg2"
adopath +  "N:/SIL-Common/reghdfe-master/package"

log using "./4.PS_leadlag_zip.log", replace

	program getColnames 
	 quietly {
	  matrix A = e(b)'
	 
	 local rownames : rowfullnames A
	 local c : word count `rownames'
	 local names : colfullnames A

	 svmat A, names(col)

	 gen rownames = ""
	 forvalues i = 1/`c' {

	  replace rownames = "`:word `i' of `rownames''" in `i' 
	   } 

	 }

	end


	program getRegStats
	 quietly {
	  matrix A = e(b)'
	 
	 local rownames : rowfullnames A
	 local c : word count `rownames'
	 local n = `c' + 1
	 local pll = `c' + 2
	 local rsq = `c' + 3
	 local ndrop = `c' + 4

	 replace y1 =e(N) in `n' 
	 replace y1 =e(ll) in `pll' 
	 replace y1 =e(r2_p) in `rsq' 
	 replace y1 =e(N_drop) in `ndrop' 

	 }

	end


	use ./MatchedSample_Full.dta 
 
 ** START STATISTICS **

*Combined sexes
 
 levelsof leadlag if leadlag!=16 & leadlag~=0, local(lllocal)

		foreach var of varlist oc_glaucoma oc_diabcomp1 oc_diabcomp2 {
		*This limits the sample to diabetes patients
		 forval j=3(1)3  {

		if regexm("`var'","oc_([a-z0-9]+)$") { 
		  local pre = "`=regexs(1)'"

		 capture drop zip3`pre' zipct`pre'
		 gen zip3`pre'=zip3dig 
		 egen zipct`pre'=count(zip3dig) if kv_`pre'==1 , by(zip3dig)
		 replace zip3`pre'=0 if zipct`pre'<1250
		 replace zip3`pre'=0 if kv_`pre'==0 
		 replace zip3`pre'=STATE_encode+1000 if zip3`pre'==0 
		 
		 capture rename mxicd3_chronic mxicd3_endocrine 

		logit `var' i.sex_pat#ib16.leadlag vmerger##c.nyear i.vmerger##i.sex_pat  nonwhite#i.sex_pat i.p_age nonwhite#i.agecat i.agecat##i.sex_pat i.qtr  ib452.zip3`pre' fmp pspcl30 pspcl06 pspcl69 pspcl18 pspcl48 pspcl93 pspclOther mxicd1_acute mxicd5_acute mxicd18_acute mxicd*_c* fs* lnvisit* [pweight=samp_weight] if flag`j'==1 & kv_`pre'==1 , cluster( propScoreId)noomit

			foreach i of local lllocal {
			 
				nlcom (exp(_cons +  _b[1.sex_pat#`i'.leadlag]+_b[1.vmerger] + _b[77.p_age])/(1+exp(_cons +_b[1.sex_pat#`i'.leadlag]+_b[1.vmerger] + _b[77.p_age])) - (exp(_cons + _b[1.vmerger] + _b[77.p_age])/(1+exp( _cons +_b[1.vmerger] + _b[77.p_age]))))
				 
				  mat A`i'=r(b)
				  mat B`i'=r(V)
				  mat C`i'=[A`i',B`i']
				  
				   if `i'==1 {
				   mat D=C`i''
				   }
				   if `i'!=1 {
				   mat D=[D,C`i'']
				   }
			  }

			foreach i of local lllocal {
				
				nlcom (exp(_cons +_b[2.sex_pat#`i'.leadlag]+_b[1.vmerger] +_b[2.sex_pat] + _b[77.p_age] + _b[3.agecat#2.sex_pat])/(1+exp(_cons +_b[2.sex_pat#`i'.leadlag]+_b[1.vmerger] +_b[2.sex_pat] + _b[77.p_age] + _b[3.agecat#2.sex_pat])) - (exp( _cons +_b[1.vmerger] +_b[2.sex_pat] + _b[77.p_age] + _b[3.agecat#2.sex_pat])/(1+exp( _cons +_b[1.vmerger] +_b[2.sex_pat] + _b[77.p_age] + _b[3.agecat#2.sex_pat])))) 
				
				  mat A`i'=r(b)
				  mat B`i'=r(V)
				  mat C`i'=[A`i',B`i']
				  if `i'==1 {
				   mat E=C`i''
				   }
				   if `i'!=1 {
				   mat E=[E,C`i'']
				   } 
			  }
			  
		  matrix beta  = e(b)' 
		  matrix sesq = vecdiag(e(V))' 
		  svmat beta
		  getColnames
		  getRegStats
		  svmat sesq
		  rename rownames rn`pre'_`j'
		  capture drop y1 
		  capture drop j 
		  capture drop i
		  rename beta1 beta_`pre'_`j'
		  rename sesq1 sesq_`pre'_`j'
		  matrix Dt = D'
		  matrix Et = E'
		  svmat Dt
		  svmat Et
		  rename Dt1 mem_beta_`pre'_`j'
		  rename Dt2 mem_sesq_`pre'_`j'
		  rename Et1 mef_beta_`pre'_`j'
		  rename Et2 mef_sesq_`pre'_`j'
		  
		   }
		  }
		}

		save "./LL_diab_zip.dta", replace

		use "./LL_diab_zip.dta" 
		rename mxicd3_endocrine mxicd3_chronic 

		levelsof leadlag if leadlag!=16 & leadlag~=0, local(lllocal)

		foreach var of varlist oc_acuteicd8 oc_ihd oc_ami {
		*This limits the sample to hyptertensive patients
		 forval j=2(1)2  {

		if regexm("`var'","oc_([a-z0-9]+)$") { 
		  local pre = "`=regexs(1)'"

		 capture drop zip3`pre' zipct`pre'
		 gen zip3`pre'=zip3dig 
		 egen zipct`pre'=count(zip3dig) if kv_`pre'==1 , by(zip3dig)
		 replace zip3`pre'=0 if zipct`pre'<1250
		 replace zip3`pre'=0 if kv_`pre'==0 
		 replace zip3`pre'=STATE_encode+1000 if zip3`pre'==0 

		 capture rename mxicd8_chronic mxicd8_circulatory
		  
		logit `var' i.sex_pat#ib16.leadlag vmerger##c.nyear i.vmerger##i.sex_pat  nonwhite#i.sex_pat i.p_age nonwhite#i.agecat i.agecat##i.sex_pat i.qtr  ib452.zip3`pre' fmp pspcl30 pspcl06 pspcl69 pspcl18 pspcl48 pspcl93 pspclOther mxicd1_acute mxicd5_acute mxicd18_acute mxicd*_c* fs* lnvisit* [pweight=samp_weight] if flag`j'==1 & kv_`pre'==1 , cluster( propScoreId)noomit
		 

			foreach i of local lllocal {
			 
				nlcom exp(_cons +  _b[1.sex_pat#`i'.leadlag]+_b[1.vmerger] + _b[77.p_age])/(1+exp(_cons +_b[1.sex_pat#`i'.leadlag]+_b[1.vmerger] + _b[77.p_age])) - (exp(_cons + _b[1.vmerger] + _b[77.p_age])/(1+exp( _cons +_b[1.vmerger] + _b[77.p_age])))
				 
				  mat A`i'=r(b)
				  mat B`i'=r(V)
				  mat C`i'=[A`i',B`i']
				  
				   if `i'==1 {
				   mat D=C`i''
				   }
				   if `i'!=1 {
				   mat D=[D,C`i'']
				   }
			  }

			foreach i of local lllocal {
				
				nlcom (exp(_cons +_b[2.sex_pat#`i'.leadlag]+_b[1.vmerger] +_b[2.sex_pat] + _b[77.p_age] + _b[3.agecat#2.sex_pat])/(1+exp(_cons +_b[2.sex_pat#`i'.leadlag]+_b[1.vmerger] +_b[2.sex_pat] + _b[77.p_age] + _b[3.agecat#2.sex_pat])) - (exp( _cons +_b[1.vmerger] +_b[2.sex_pat] + _b[77.p_age] + _b[3.agecat#2.sex_pat])/(1+exp( _cons +_b[1.vmerger] +_b[2.sex_pat] + _b[77.p_age] + _b[3.agecat#2.sex_pat])))) 
				 
				  mat A`i'=r(b)
				  mat B`i'=r(V)
				  mat C`i'=[A`i',B`i']
				  if `i'==1 {
				   mat E=C`i''
				   }
				   if `i'!=1 {
				   mat E=[E,C`i'']
				   } 
			  }

			  
		  matrix beta  = e(b)' 
		  matrix sesq = vecdiag(e(V))' 
		  svmat beta
		  getColnames
		  getRegStats
		  svmat sesq
		  rename rownames rn`pre'_`j'
		  capture drop y1 
		  capture drop j 
		  capture drop i
		  rename beta1 beta_`pre'_`j'
		  rename sesq1 sesq_`pre'_`j'
		  matrix Dt = D'
		  matrix Et = E'
		  svmat Dt
		  svmat Et
		  rename Dt1 mem_beta_`pre'_`j'
		  rename Dt2 mem_sesq_`pre'_`j'
		  rename Et1 mef_beta_`pre'_`j'
		  rename Et2 mef_sesq_`pre'_`j'
		  

		   }
		  }
		}

		save "./LL_hyp_zip.dta", replace

		foreach var of varlist oc_death {
		 forval j=1(1)1  {


		 if regexm("`var'","oc_([a-z0-9]+)$") { 
		  local pre = "`=regexs(1)'"

		 capture drop zip3`pre' zipct`pre'
		 gen zip3`pre'=zip3dig 
		 egen zipct`pre'=count(zip3dig) if kv_`pre'==1 , by(zip3dig)
		 replace zip3`pre'=0 if zipct`pre'<1250
		 replace zip3`pre'=0 if kv_`pre'==0 
		 replace zip3`pre'=STATE_encode+1000 if zip3`pre'==0 

		 capture rename mxicd8_chronic mxicd8_circulatory
		  
		logit `var' i.sex_pat#ib16.leadlag vmerger##c.nyear i.vmerger##i.sex_pat  nonwhite#i.sex_pat i.p_age nonwhite#i.agecat i.agecat##i.sex_pat i.qtr ib452.zip3`pre' fmp pspcl30 pspcl06 pspcl69 pspcl18 pspcl48 pspcl93 pspclOther mxicd1_acute mxicd5_acute mxicd18_acute mxicd*_c* fs* lnvisit* [pweight=samp_weight] if flag`j'==1 & kv_`pre'==1 , cluster( propScoreId)noomit
		 

			foreach i of local lllocal {
			 
				nlcom exp(_cons +  _b[1.sex_pat#`i'.leadlag]+_b[1.vmerger] + _b[77.p_age])/(1+exp(_cons +_b[1.sex_pat#`i'.leadlag]+_b[1.vmerger] + _b[77.p_age])) - (exp(_cons + _b[1.vmerger] + _b[77.p_age])/(1+exp( _cons +_b[1.vmerger] + _b[77.p_age])))
				 
				  mat A`i'=r(b)
				  mat B`i'=r(V)
				  mat C`i'=[A`i',B`i']
				  
				   if `i'==1 {
				   mat D=C`i''
				   }
				   if `i'!=1 {
				   mat D=[D,C`i'']
				   }
			  }

			foreach i of local lllocal {
				
				nlcom exp(_cons +_b[2.sex_pat#`i'.leadlag]+_b[1.vmerger] +_b[2.sex_pat] + _b[77.p_age] + _b[3.agecat#2.sex_pat])/(1+exp(_cons +_b[2.sex_pat#`i'.leadlag]+_b[1.vmerger] +_b[2.sex_pat] + _b[77.p_age] + _b[3.agecat#2.sex_pat])) - (exp( _cons +_b[1.vmerger] +_b[2.sex_pat] + _b[77.p_age] + _b[3.agecat#2.sex_pat])/(1+exp( _cons +_b[1.vmerger] +_b[2.sex_pat] + _b[77.p_age] + _b[3.agecat#2.sex_pat])))
				 
				  mat A`i'=r(b)
				  mat B`i'=r(V)
				  mat C`i'=[A`i',B`i']
				  if `i'==1 {
				   mat E=C`i''
				   }
				   if `i'!=1 {
				   mat E=[E,C`i'']
				   } 
			  }

			  
		  matrix beta  = e(b)' 
		  matrix sesq = vecdiag(e(V))' 
		  svmat beta
		  getColnames
		  getRegStats
		  svmat sesq
		  rename rownames rn`pre'_`j'
		  capture drop y1 
		  capture drop j 
		  capture drop i
		  rename beta1 beta_`pre'_`j'
		  rename sesq1 sesq_`pre'_`j'
		  matrix Dt = D'
		  matrix Et = E'
		  svmat Dt
		  svmat Et
		  rename Dt1 mem_beta_`pre'_`j'
		  rename Dt2 mem_sesq_`pre'_`j'
		  rename Et1 mef_beta_`pre'_`j'
		  rename Et2 mef_sesq_`pre'_`j'
		  

		   }
		  }
		}
		 
		 
		save "./LL_death_zip.dta", replace
		clear
		 